#!/usr/bin/python

import re
import gc
import sys
import bsddb
import string
import itertools

BORDERS_RE1 = re.compile('(complement)?\(?<?(\d+)..>?(\d+)\)?')
BORDERS_RE2 = re.compile('(complement)?\(?join\((.+),(.+)\)\)?')

TABLE = string.maketrans('actg', 'tgac')

def complement(s):
    return s.translate(TABLE)

def reverse_str(s):
    return ''.join(reversed(s))

def transparent(s):
    return reverse_str(complement(s))

def make_gene_simple(genom, simple_borders):
    m = BORDERS_RE1.match(simple_borders)
    if m:
        is_compl, left, right = m.groups()
        is_compl = bool(is_compl)
        left     = int(left)
        right    = int(right)
        gene = genom[left-1:right]
        if is_compl:
            gene = transparent(gene)
        return gene
    raise Exception('Bad simple borders: %s' % `simple_borders`)

def make_gene(genom, borders):
    m = BORDERS_RE1.match(borders)
    if m:
        return make_gene_simple(genom, borders)
    
    m = BORDERS_RE2.match(borders)
    if m:
        is_compl, borders1, borders2 = m.groups()
        gene = make_gene_simple(genom, borders1) + make_gene_simple(genom, borders2)
        if is_compl:
            gene = transparent(gene)
        return gene
    

if __name__ == "__main__":
    gc.set_threshold(700, 1000, 10) 

    genom_db_fname = sys.argv[1]
    gen_desc_db_fname = sys.argv[2]

    cnt_all = 0
    cnt_genes = 0
    genom_db = bsddb.hashopen(genom_db_fname, 'r')
    gen_desc_db = bsddb.hashopen(gen_desc_db_fname, 'r')
    
    for locus, grp in itertools.groupby(sorted(gen_id.rsplit('_', 1) for gen_id in gen_desc_db.iterkeys()), key=lambda i: i[0]):
        cnt_all += 1
        genom = genom_db[locus]
        for locus, gene_md5 in sorted(grp):
            gene_id = "%s_%s" % (locus, gene_md5)
            borders = gen_desc_db[gene_id].split('\n')[0]
            gene = make_gene(genom, borders)
            if gene:
                print "%s\t%s" % (gene_id, gene)
                cnt_genes += 1
            else:
                print >> sys.stderr, 'Unknown borders: %s %s' % (gene_id, `borders`)
    genom_db.close()
    gen_desc_db.close()
    gen_desc_db.close()

    print >> sys.stderr, "Done! (%s/%s)" % (cnt_all, cnt_genes)

